%% Simulates set of basic shocks
% Read the model
[m,p,mss] = readmodel(false);

%% Shocks
% List of shocks
% One period unexpected shocks: inflation, output, exchange rate, interest rate
% Creates a list (vector) of shocks and list of their names. Used shocks' names must
% be the names found in the model code (in file 'model.model') except of the prefix
% 'e_'
listshocks = {'e_dot_cpi_x','e_dot_wfood','e_dot_woil','e_lgdp_gap','e_ls','e_rn'};
listtitles = {'Inflationary Shock','World Food Prices Shock','Wordl Oil Prices Shock',...
    'Aggregate Demand Shock', 'Exchange Rate Shock', 'Interest Rate Shock'};

% Sets the time frame for the simulation 
startsim = qq(1,1);
endsim = qq(4,4); % four-year simulation horizon

% For each shock a zero database is created (command 'zerodb') and named as 
% database 'd.{shock_name}'
for i = 1:length(listshocks);
    d.(listshocks{i}) = zerodb(m,startsim-4:endsim);
end

% Fills the respective databases with the shocks' values for the starting
% point of the simulation (startsim). For simplicity, all shocks are set to
% 1 percent
d.e_dot_cpi_x.e_dot_cpi_x(startsim) = 1;
d.e_dot_wfood.e_dot_wfood(startsim) = 1;
d.e_dot_woil.e_dot_woil(startsim) = 1;
d.e_lgdp_gap.e_lgdp_gap(startsim) = 1;
d.e_ls.e_ls(startsim) = 1;
d.e_rn.e_rn(startsim) = 1;

% Simulates the shocks using the command 'simulate'. For this model
% 'm' and respective database 'd.{shock_name} are used. Results are written 
% in object 's.{shock_name}'. Command 'dboverlay' merges the historical 
% database d.{shock_name} with simulation database s.{shock_name}. 
p = plan(m,startsim:endsim);
for i=1:length(listshocks);    
    s.(listshocks{i}) = simulate(m,d.(listshocks{i}),startsim:endsim,'deviation',true);
    s.(listshocks{i}) = dbextend(d.(listshocks{i}),s.(listshocks{i}));
end

% Report
x = report.new('Shocks','visible',true);

% Figures
sty = struct();
sty.line.linewidth = 1.5;
sty.line.linestyle = {'-';'--'};
sty.line.color = {'k';'k'};
sty.axes.box = 'on';
sty.legend.location = 'southOutside';
sty.legend.orientation = 'horizontal';

for i = 1:length(listshocks);
x.figure(listtitles{i},'subplot',[3,3],'style',sty,'range',startsim:endsim,'dateformat','YYYY:P');

x.graph('Inflation','legend',true);
x.series('',[s.(listshocks{i}).dot_cpi s.(listshocks{i}).dot4_cpi],'legendEntry=',{'q-o-q','y-o-y'});

x.graph('Core Inflation','legend',false);
x.series('',[s.(listshocks{i}).dot_cpi_x]);

x.graph('Food Inflation','legend',false);
x.series('',[s.(listshocks{i}).dot_cpi_food]);

x.graph('Oil Inflation','legend',false);
x.series('',[s.(listshocks{i}).dot_cpi_oil]);

x.graph('Nominal Interest Rate','legend',false);
x.series('',s.(listshocks{i}).rn);

x.graph('Nominal Exchange Rate','legend',true);
x.series('',[s.(listshocks{i}).dot_s (s.(listshocks{i}).ls - s.(listshocks{i}).ls{-4})],'legendEntry=',{'q-o-q','y-o-y'});

x.graph('Output Gap','legend',false);
x.series('',[s.(listshocks{i}).lgdp_gap]);

x.graph('Monetary Conditions','legend',false);
x.series('', [s.(listshocks{i}).mci]);

x.graph('RIR and RER Gaps','legend',true);
x.series('', [s.(listshocks{i}).rr_gap s.(listshocks{i}).lz_gap],'legendEntry=',{'RIR', 'RER'});

x.pagebreak();

end

x.publish('Shocks','display',false);